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Abstract 

We propose a hybrid model for optimal control of a bioreactor which features a full coupling 
between a linear optimization problem and a dynamic one. 

1 Introduction 

Biofuels provide a concrete answer to the pressing need for renewable energy. The problem of 
increasing the efficiency and reducing the cost of biofuel production has been subject of intensive 
research [2l|9]. A feasible source to obtain biofuels consists in using ethanol produced by cyanobac- 
teria and microalgae [121 US] • This paper deals with the application of optimal control techniques 
to a general bioreactor for biofuel production. 

Since the last decades, optimization of bioprocesses has been a line of research connecting opti- 
mal control to system biology, see [24j and references therein. Progress in plant genetic engineering 
has opened novel opportunities to use plants as bioreactors for safe and cost effective production 
of vaccine antigens. A review of methods and applications of plant, tissue and cell culture based 
expression strategies and their use as bioreactors for large scale production of pharmaceutically im- 
portant proteins can be found in [23j. Exploiting bioreactors was also proved a fruitful method to 
deal with the problem of wastewater treatment in environmental engineering [19] . In this context a 
dynamic optimization problem arising frequently is the minimization of the time needed to reach a 
fixed target configuration for the bioreactor, see [22J. In this case, the trajectory evolves according 
to a system of nonlinear ODEs and may satisfy some boundary constraints as well. The right-hand 
side of the evolution equations involves not only control functions, i.e. parameters through which 
we can modify the dynamics, but also some unknown functions (e.g. biomass growth rate) of the 
evolving quantities themselves. Therefore, computing optimal controls through standard techniques 
(necessary conditions) is affected by the unknown functions' behavior. An interesting approach to 
deal with this issue using observability techniques was proposed in f31, where an application to 
bioreactors is also provided (see also [llj). 

Other applications of optimal control techniques have been studied in the field of biological 
systems [Ij. In the context of medical treatment, chemotherapy was used |16j as a dynamic control 
to optimize treatment scheduling in early stage HIV-infected cases. In that paper the authors 
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use the effect of chemotherapy on viral production to maximize benefits in terms of T cell count 
and minimize the systemic cost of the treatement. Time dependent control strategies were also 
exploited in models to contain the emergence of drug-resistant strains of tubercolosis [15]. Here 
controls are represented by efforts in finding patients in which virus is only latent and in completing 
treatment for patients in which virus is already active. The objective function balances the effect 
of minimizing the cases of latent and infectious drug-resistant tubercolosis and minimizing the cost 
of implementing the control treatments. 

In this paper we are concerned with the optimization of a bioprocess for biofuel production. 
Namely, we consider a model for the metabolic activity of a microorganism (e.g. E. Coli or Saccha- 
romyces Cerevisiae) in a fed-batch culture with different feeding substrates. The optimal control 
problem is to maximize the productivity of a certain side metabolite (e.g., ethanol) through different 
feeding rates. 

A steady-state approach to model cellular metabolism is Flux Balance Analysis (FBA, see [20j). 
The main assumption of FBA is that metabolic activity of cells is performed in such a way that 
the growth rate of cells is maximized. Since genome-scale stoichiometric models for bacteria such 
as E. Coli are available, this translates in a linear program where the objective is cellular growth 
rate and constraints are given by metabolic reactions. Optimization is then performed by means 
of genetic manipulations on the bacteria (gene deletions and insertions). Other approaches based 
on flux balance analysis were proposed that take account of transcriptional and regulatory effects 
in [7|, that couple the steady-state metabolic activity with a dynamic model [171 [T8] and that 
integrate both aspects [8j. Optimizing ethanol productivity in fed-batch cultures modeled through 
Dynamic Flux Balance Analysis [H] consists in using outputs of FBA (cellular growth rate and 
metabolite fluxes) to update at each time step the dynamics for evolving extracellular quantities. 
Within this framework, control can be performed at two levels: intracellular controls (genetic 
modifications), which are implemented by acting on constraints of the FBA, and extracelullar 
controls (of dynamic nature), which are implemented by changing feeding rates. Genetic strategies 
were deeply investigated in [14] for in silico evolution of a yeast strain in glucose and xylose media 
to maximize ethanol productivity with constant feeding rates. 

In this paper, we analyze extracellular controls, i.e., feeding rates for glucose and xylose, which 
we allow to be time dependent. Namely, we cast the problem of optimizing ethanol for the biore- 
actor in the language of optimal control, where concentrations of feeding substrates are modified 
by an external agent in order to control ethanol production. Even though we do not consider 
genetic manipulations, our model can integrate intracellular controls. Our purpose is to develop 
up a general model which shares a full coupling between extracellular dynamics and intracellular 
metabolic activity: at each time instant solutions of ODEs feed fluxes to constraint the FBA, 
whereas FBA provides cellular growth rate and ethanol uptake ODEs see Figure [TJ Since both 
the objective function and the constraints in FBA are linear (with respect to fluxes), outputs are 
piecewise linear. Dependence of fluxes on metabolite concentrations is modeled through Michaelis- 
Menten behavior, i.e., through rational functions. Therefore, FBA outputs appearing in dynamics 
are piecewise smooth functions of evolving quantities. The main idea is to implement the coupling 
between the intracellular and extracellular level through a hybrid control system where in each 
location the outputs of FBA are smooth. One of the advantages of a hybrid formulation is that it 
allows to take account of different timescales. One can either assume that information translates 
instantaneously from the micro level to the macro one (and viceversa) or one can implement delays, 
which are observed in experimental data. This is useful as there are configurations of the extra- 
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Figure 1: Bioprocess scheme exhibiting full coupling between metabolic activity and external dy- 
namics 



cellular environment (such as saturation of glucose or oxygen) which may cause a major change in 
the metabolic pathway involved inside the cell with a consequent delay in the variation of outputs 
of the FBA. 

As a first step toward a unified model, we focus on the analysis of optimal trajectories contained 
in a connected component where FBA parameters are smooth. In other words, we study what 
happens in each location of the hybrid model. More precisely, the location is characterized by a 
control-affine system of the type 

X = Fo{x,y'') + uiFi{x) + U2F2{x), 

where a; is a vector containing substrate concentrations (biomass, glucose, xylose, ethanol), ui,U2 
are feed rates of glucose and xylose, are parameters coming from FBA, and F(),Fi,F2 are 
smooth. We consider an optimal control problem in Mayer form where we optimize the total 
amount of ethanol at the final time. Since the system is control-affine, the maximization condition 
of the Pontryagin Maximum Principle allows to compute the control along an extremal trajectory 
as a feedback law, as long as the corresponding switching function (derivative of Hamiltonian with 
respect to control) has only isolated zeros. Otherwise, when the trajectory shows a singular arc, 
i.e., a time interval where a switching function vanishes identically, the main tool to find controls is 
exploiting conditions given by annihilation of higher order derivatives with respect to time. As our 
interest is driven by applications, one would rather avoid singular arcs, which represent an obstacle 
to efficient and reliable numerical simulations. The main result of the paper describes properties of 
extremal trajectories having singular arcs, under the assumption that parameters are affine with 
respect to substrates' fluxes and that growth is anaerobic. Notice that, due to specific properties 
of the optimal control problem considered here, one cannot expect that optimal singular arcs do 
not occur, even generically. Indeed, as the vector fields Fi,F2 are actually constant, we cannot use 
results on codimension of singular trajectories, such as [6], where the authors deal with control- 
affine systems satisfying generic assumptions. We also consider the role of oxygen as a control. 
In the case where feeding rates of glucose and xylose are piecewise constant, we characterize the 
solution of the adjoint system along singular arcs. 

The structure of the paper is the following. In Section [2] we recall the general model of a 
bioreactor and formulate the optimal control problem where the cost to be maximized is ethanol 
productivity of the bioreactor (i.e., a certain function of the total amount of ethanol produced 
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at the final time.) In Section [3] we state and prove the main result concerning optimal singular 
trajectories. In Section Owe discuss the role of oxygen as a control. 



2 Problem formulation 

Bioreactors are processes where a living microorganism metabolizes some substrates and conse- 
quently grows and produces other metabolites. In [Tl| the authors consider in silico evolution of a 
yeast strain, Saccharomyces cerevisiae, which grows in a fed-batch culture with glucose and xylose 
and produces ethanol. The dynamic model can be applied to a general bioreactor. 

We denote by V the total culture volume, which is assumed to grow linearly with respect to time 
with constant rate F. Biomass concentration is denoted by X. The total biomass in the culture 
evolves linearly with a growth rate fi. Concentrations of feeding substrates (glucose G and xylose 
Z) are characterized by an evolution which takes account of a feeding rate and a compensation 
term due to metabolism of the microorganism. Feeding rates for substrates represent the control 
we perform on the system and are denoted by ui,U2- The organism metabolizes glucose and xylose 
with specific rates (or fluxes) Vg,Vz- Oxygen concentration O is also involved in the bioprocess and 
represent a further control whose role we analyze in Section HI The produced metabolite we study 
is ethanol. Denoting ethanol concentration by E, we assume that the organism produces ethanol 
proportionally to the total biomass through a specific rate Ve- Therefore, the control system we 
analyze is 

' y = F 

V'X = fi{G,Z,E,0)VX 

VG = Fui-Vg{G,E)VX (1) 

VZ = Fu2-v,{G,Z,E)VX 

VE = Ve{G,Z,E,0)VX. 

It is evident from system ([1]) the parameters fj,,Vg,Vz,Ve depend on evolving variables. First of all, 
oxygen uptake kinetics follows Michaelis-Menten kinetics 

Vo{0) ^ 



o max , I ,o ' 



whereas glucose uptake kinetics has an additional regulatory term to capture growth rate suppres- 
sion due to high ethanol concentration, i.e., 

Vg{G,E) = V, ^ ^ 



-'kg+Gl+E/kl 



Xylose uptake kinetics has a similar form with another regulatory term to account for inhibited 
xylose metabolism in presence of the preferred substrate (glucose) , 

Vz{G,Z,E)=v ^ ^ ^ 



'"^^""kz + Zl + E/kf^l + G/k 



ig 



Parameters Uomaxj 'Wgmaxj 'Wzmax, kg, kf^, kz, kf^i kig are positive and constant. As for /i, Ve, the model 
is based on the principle that the metabolic activity of the microorganism is performed so that 
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biomass growth rate is maximized. This is equivalent to say that and Ve are outputs of an 
optimization problem stated as 



n 



fi{G, Z, E, O) = max WjVj (2) 
i=i 
s.t. Sv = 

Vg = Vg{G,E) 

v, = v,{G,Z,E) 
v, = v,{G,Z,E) 

Vo = VoiO) 

<Vj <Vj, j = l,...,n 

VeiG, Z, E, O) = argmax Z, E, O). (3) 

In (l2|), -u G R"' is the vector of fluxes (among which glucose, xylose, oxygen and ethanol fluxes) 
considered in the model; w G [0, 1]" is the vector of weights which determines fluxes producing 
biomass; S E M^'^"(M) is the stoichiometric matrix which encode the metabolic network inside 
the cell (involving r reactions and n metabolites); v is an upper bound associated with the mi- 
croorganism, (see [20j for an exhaustive treatment of metabolic networks in systems biology.) For 
our purposes, w, S, v are considered as given parameters. They depend on the specific strain of 
microorganism used in the bioreactor and are usually determined through experimental data see 
for instance supplementary data of [8] for E. Coli or [13j for Saccharomyces Cerevisiae. 

Equation ([3]) is to be read as follows. Let V^^^ denote the set of vectors v satisfying the 
constraints above and realizing the maximum, i.e., n{G,Z,E,0) = "^jWjVj. Then Ve{G, Z, E,0) 
is the maximum of Vg as v varies in y™^^. We fix the final time tf > and an initial condition 
(Vo) -^0; Co, .^0) Oo) and we seek to optimize ethanol productivity along a trajectory of ([1]), i.e., 

max —. ' (4) 

(«i,«2)GW j^f{ui{s) + U2{s)) ds 

where the class of admissible controls is 

U = {{ui,U2) € L°°([0,t/])2 I Ui{t) G [0, 1] for almost every t}. 

We may assume that the functional does not diverge, as the trivial control strategy where there 
is no feed, i.e., {ui{t),U2{t)) = (0,0), is not optimal. 

The problem shares a full coupling between a classical optimal control problem ([I]), ([!]), which 
describes the dynamics outside the cell, and a linear optimization problem ([2]), which models the 
metabolic activity inside the cell. In other words, at each time instant t, on the one hand, one 
needs to solve the LP ([2]) to obtain /U,t;e at time t to plug into ([T]); on the other hand, dynamics 
in ([1]) must be integrated to in order to get the triple {G{t), Z{t), E{t)) that allows to determine 
constraints in ([2]). 

Since in ([2]) both the objective function and the constraints are linear with respect to v, the 
outputs fj, and Ve are piecewise linear with respect to v. Therefore fi, are piecewise smooth as 
functions of {Z, G, E, O). In the sequel we assume that the outputs of ([2]) are smooth. This amounts 
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to say that the system evolves in a domain where /U,Ue are smooth as functions of {Z,G,E,0). 
Under this assumption, we study optimal trajectories for problem ([I]), (jlj). 

Let us give a more compact formulation of the optimization problem above. We will deal first 
with the case where growth is anaerobic, that is, oxygen concentration O is zero (see Section |4] for 
the analysis where oxygen is taken as a control.) Let us rename 

X3 = VG, X4 = VZ, X5 = VE, X6 = VX, xj = V. (5) 

It is natural to add two new variables to ([1]) by setting 

Xi{t) = / Ui{s)ds, i = 1,2, 
Jo 

so that the control system ([T]) reads 

Xl = Ui 
±2 = U2 

X3 = Fui - Vg{x3,X5, X7)xq 

Xi = Fu2-Vz{x3,X4„X5,X7)xq (6) 
X5 = Ve{x3,X4,X5,X7)x6 
Xe = fi{x3,Xi,X5,X7)xQ 

X7 = F 
where 

{kgX7 + X3){X7 + X^/kf^) 
2 

Vz{x3,X4,X5,X7) = t'^maxT] : w , ^^^I, w ] TTTT' (8) 

{kzX7 + X4){X7 + X3/kig)[X7 + Xs/fcfJ 

and t'e(x3, X4, X5, X7), ^(2:3, X4, X5, X7) are defined in the obvious way by (l2|), ([3]) using ([5]). 

The equation for xy can be trivially integrated. Indeed, the total volume of the culture plays 
the same role as time along the experiment. Here we prefer to keep the state variable X7 in order 
the control system to be autonomous. 

Setting X = (xi, X2, X3, X4, X5, xg, X7), the optimization in (jH) becomes 

max (9) 

{U1,U2)&A 

where ip{x) = ^^^^^ ■ The latter is smooth on the subset 7^ = {x € | xi + X2 > 0}. Every 
trajectory starting at a point with xi(0) + X2(0) > satisfies the state constraint xi{t) + X2{t) > 0, 
i.e., it belongs to the region TZ for every t > 0. 

The control system ([6]) is affine with respect to controls and can be written in the form 

X = /(x, n), X € M"^ 
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where u = (ni,U2), f{x,u) = Fq{x) + uiFi{x) + U2F2{x) and 



Foix) 






-Vg{x3,X5,X7)xG 
-Vzix3,X4,X5,X7)xQ 
Ve{x3,X4,X5,X7)x6 
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We end this section by stating a result which ensures the existence of optimal solutions to ([T]) , 
with suitable initial conditions. 



Proposition 1 Assume v^, n G C°°(M'',M+) and let x^ belong to 

V = {x \ Xi>0,i = 1. . .6,xi + X2 > 0,X7 > 0}. 
Then there exist an optimal solution to ^ satisfying the initial constraint x(0) 



(10) 



Sketch of proof. Under the assumption that fi, Ve are smooth, every trajectory of the control 
system starting at a point x^ in the set (jlOp . is well-defined for all t > and satisfies Xi{t) > for 
every t > 0. Indeed, let T>q = {x G M7 \ Xi > 0,i = 1 . . . Q, xj > Xj}. Then the dynamics is smooth 
and satisfies the sublinear bound 



|/(x,n)| <C(l + |x|), VxGPo, (ni,n2) G [0,1]2, 

for a certain constant C > 0. Let T be such that a trajectory starting at x^ is well-defined on 
[0,T]. Then x{t) G T>q for every t G [0,T]. Indeed, inequalities Xi{t) > for i = 1,2, and 
xjit) > are trivially satisfied. Similarly, since the equation for xq is linear with respect to xq, 
we have XQ{t) > 0. By assumption (and thanks to the constraints < Vj < vj in ([2])), ethanol fiux 
'ye(x3, X4, X5, xy) is non negative for every t, whence X5(t) > 0. Finally, let t be such that X3(t) = 0. 
Then X3{t) = Fu2{t) > 0, hence x^lt) > for t > t. In the same way one infers that X/^{t) > for 
every t. Hence trajectories starting at a point of the set (jlOp are well-defined for every t > and 
x(t) G Vq for all t > 0. Moreover, because x^ + X2 > 0, we have xi(tf) + X2{tf) > 0, that is, the 
final point of any trajectory belongs to the TZ on which the cost functional @ is smooth. 

Since the dynamics ([6]) is control affine, the set of velocities {f{x,u) | n G [0, 1]^} is convex for 
every x. Therefore, the existence of an optimal solution is a consequence of classical results (see [3l 
Theorem 5.1.1]). ■ 

3 Analysis of extremal trajectories 

In this section we study optimal trajectories of (l6|), ([9]) starting from a given initial condition 
x^ GV. The Hamiltonian associated with the optimal control problem is 

H{x,X,u) = (A,Fo(x))+ni (A,Fi(x))+U2 (A,F2(x)), 

where A = (Ai, A2, A3, A4, A5, Ag, A7) G MJ denotes the covector. 
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Let u*{t) = {ul{t) , U2{t)) be an admissible control such that the corresponding trajectory of 
([6]) starting at x^, denoted by x{t), is optimal. Applying the Pontryagin Maximum Principle j21j . 
there exists a solution A(t) 7^ of the adjoint system 



dH ^ 
dx 



x{t),X{t),u*it)), A(t/) = VV(x(V)), 



which reads 



dVg 

8x3 



^6 { ast A3 + it A4 



9x5 



Ai = 
A2 = 

A3 

A4 

A5 

Ae = VgXs + VzXa: 
such that, for almost every t £ [0,tf] 



dVe 



8x3 ^ 
8fj. 



X6[q±X3 + 



dvz \ , 



8Ve 



X. 



VeX^ - fiXe 



Xo -L \^ 



&7A6 



dfj. \ 
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H{xit),Xit),u*{t)) 



max (A(t), Fo{x{t))) + toi (A(t), Fi{x{t))) + (A(t), F2{x{t))) 
< cji < 1 

< LJ2 < 1 



(see also [3l Theorem 6.1.1]). The transversality condition for the covector at final time is 



X{tf) 



0,0, 



(xi(t/)+X2(t/))2' (X.itf) + X2itf)y 

and, since the Hamiltonian does not depend on xi,X2, we obtain 



Xx{tf)^X2{tf) 



Ai(t) = A2(t) 
Define the switching functions 



(^l(t/)+X2(t/)) 



2 • 



(12) 



^i{t) = {X{t),F,{x[t))) , i = l,2. 

Since the Hamiltonian is affine with respect to controls, for every t such that (pi{t)(p2{t) 7^ the 
control u*{t) is uniquely determined by u*{t) = i±£l^^iW^ i = 1,2. 

We say that an extremal x{-) has a singular arc on [a,b] C [0, t/], with a < 6, if ^i\[a^b] = 
or (/?2|[a,6] = 0. In general more degenerate singular arcs can occur: for instance, arcs where one 
among the switching functions vanishes identically and the other one has non-isolated zeros and 
does not vanish identically. In this case, the corresponding component of the control features an 
accumulation of switchings. Accumulations of zeros of the switching functions can indeed happen 
along an optimal trajectory and they are known as chattering or Fuller's phenomena (see [10].) 
In this paper we assume that chattering phenomena do not occur. Even though such pathological 
phenomenon can happen, novel techniques have been recently developed to overcome it, see [5]. 

When seeking an optimal synthesis, a fundamental step is to reduce the problem to a finite- 
dimensional one. For instance, this can be done if, by means of necessary conditions, one deduces 
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that every extremal trajectory is a finite concatenation of bang arcs, where the control has a simple 
form. In general this is hard to prove, still one can expect that only "a few" extremal trajectories 
are "bad", i.e., they do not satisfy this property. To do so, a first step is the result below, which 
describes the covector along badly behaved extremal trajectories for the optimal control problem 

Concerning the presence of optimal singular trajectories for control-affine systems, a genericity 
result has been shown in [6]. In that paper, the authors prove that, given a quadratic cost, generi- 
cally with respect to the control system there do not exist nontrivial optimal singular trajectories 
[6l Corollary 2.9]. Notice that this result does not apply here because system ^ is not generic, 
Fi,F2 being constant vector fields. 

Theorem 1 Assume that fj,,Ve are affine with respect to {vg,Vz), i.e., 

H = aiVg + a2Vz + fl 

Ve = hVg + b2Vz + V, (13) 

with oi, 02, 6i, &2i A) positive constant such that 

aib2 - 0261 / 
ih - b2)fi + ia2 - ai)v / 
bip, — aiv ^ 0. 

Let x{-) be an optimal trajectory having a singular arc on [a,b] C [0,tf] such that the sets of zeros 
of Lfi and ip2 are finite unions of intervals. Then one among ipi,(p2 does not vanish identically on 
[a,b]. 

If a singular arc of type = 0, (/32 7^ is optimal, then U2 = and 

• if (A, [Fi, [Fq, Fi]]) = then A4(t) is uniquely determined by x{-) and x{tf) and Xj is constant 
for every j 7^ 4. 

• if (A, [-Fi, [Fo,-Fi]]) has only isolated zeros then is uniquely determined as a feedback. 

Proof. Since there is no chattering, on [a, b] either both switching functions vanish identically or 
one is identically zero and the other one has only isolated zeros. 

Assume first that 921 = and 932|[a,fe] = 0- Then 

Xiit) + FX3{t) = 
^2{t) + FX4{t) = 0. 

Jointly with p2]) . we deduce that 

A3(t) - Ht)^-^. 

Therefore, Ai(t), A2(t), A3(t), A4(t) are uniquely determined as functions of x{tf). In the following 
we set Xi{t) = Xi{tf),X2{t) = Xi{tf),X3{t) = -Ai(t/)/F, A4(t) = -Xi{tf)/F. 
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Conditions ipi = 0,ip2 = imply jointly with the third and fourth equation in (jlip that the pair 
(A5, Ag) must satisfy 

yuTi V'^ey y — F^gsi / 

System (fT4]l can be seen as a linear system in {\^{t) , \Q{t)) whose coefficient depend only on x{t) 
and on x{tf). When fi,Ve are given as in (jl3p . the determinant of ()14p is a positive multiple of 
0162 — 02^1 • Hence, by assumption, it never vanishes and (jl4p admits a unique solution. In other 
words, we can express (X^it), XQ{t)) as a functions of x{t) and x{tf ). In what follows, we denote by 
X5(t),Xe{t) the solution of Consider now the system 

y5i = {X,[Fo,[Fo,Fi]])+ui{X,[Fi,[Fo,Fi]])+U2{X,[F2,[Fo,F,]])=0, (15) 
¥32 = (A,[Fo,[Fo,F2]])+ui(A,[Fi,[Fo,F2]])+n2(A,[F2,[Fo,F2]]) = 0. (16) 



Taking A = (Ai, A2, A3, A4, A5, Ae, A7), due to the form of fj,, Ve a computation shows that 

{X,[Fi,[Fo,Fi]])=0, (A,[F2,[Fo,Fi]]) =0, (A, [Fi, [Fq, F2]]) = 0, {X,[F2,[Fo, F^]]) = 0. 
Hence, in order system (llSp to be compatible, one needs 

(A, [Fo, [Fo, Fi]]) = 0, (A, [Fo, [Fq, F2]]) = 0. 

Nevertheless, we have 

a2Fkf^kigkzVzms.^XQX^{biX4, - 62 Aa)^ + (a2A3 - aiX4)v) 



(A,[Fo,[Fo,F2]]) 



(0261 - aib2){x5 + kf^X7){x3 + kigx-j){xi + kzX-jf 



Hence, since A3 = A4 7^ and 0162 — a2&i 7^ 0, (61 — 62)/^ + (02 — ai)v 7^ 0, we get a contradiction. 
Therefore, there cannot exist optimal singular arcs where both ipi and '^2 vanish identically. 

Assume now that <^i|[a,6] = and 932 (i) has only isolated zeros in [0,6]. In particular, 

Ut) - > 0, 

and U2{t) = i±£lli^^2(t)2 ^ thai is, U2 is given as a feedback law. Let Xi{t) = Xi{t f) ,X-i{t) = 
— Ai(t/)/F. From condition ipi = 0, which is equivalent to A3 = 0, we obtain that 



kk^,, Xs-biX^-aiXe -A4 + 62 A5 + a2Xe = n(^7) 

Setting 



A - J, 1,9., A3 - biAs - aiAs B-k'k- V -A4 + h2X^ + a2A6 

- ^9he^9r... ^ k^X^fix, + fcf.X^) ' " ^-^^^^^--(^g + fc^^^,)(^3 + k,^X,Y{x, + fe.Xr) ' 

equation (fT7|l becomes 

A + X45 = 0. 
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Let c(t) = (A, [-Fi, [-Foi^i]]) be the coefficient of ui in condition ipi = 0, see (|15p . A computation 
sfiows tfiat 

c{t) = -2F^X6X? + X, ) . (18) 

Since tliere is no cfiattering, zeros of c(-) are finite union of intervals. Hence, we distinguish two 
cases: first c\^af>\ = 0) second c(-) has only isolated zeros in [a,b]. In the second case condition 
(fi=0 allows to determine ul{t), for all t such that c{t) ^ 0, by 

, _ (A, [F2,[Fo,Fi]]) + (A,[Fo,[Fo,Fi]]) 
(A,[Fi,[Fo,Fi]]) 

Let now c(t) = 0. In this case, from ()17p and (llSp and since A;^ 7^ /cjc,, we deduce that A = and 
3:45 = 0, which are equivalent to 

fA3-6iA5(t)-aiA6(t) = 
\x4(t)(-A4(i) + 62A5(i) + a2A6(t)) = 

Using these relation and the behavior of /U,r;e, the dynamics for A5 in (jlip gives A5 = 0. Hence, 
differentiating the first equation with respect to time we obtain that Ag = 0. Therefore, imposing 
that the right hand side in the dynamics for Ag vanishes identically we obtain that the pair (A5, Ag) 
is constant. Indeed, it is the unique solution to the linear system 




d 



A = 0, 



whose coefficients depending only on x, A3 and which has determinant a\ (\)\[i—a\v). By assumption, 
this determinant is nonzero and the unique solution of the system is 

A5(t)^A3 ,_^ _ , A6(t)^-A3 



b\\x — a\v bifi — aiv 

Assuming that the zeros of X4 and B are finite union of intervals, then [a,b] = UiLi[ai,bi], where 

ai+i = bi and on each [a,, bi], Xi\\a^^bi] = or B\]^a„hi] = 0. 

Let a;4|[a._;,.] = 0. Then plugging A5, Ag into (fTTj) . the equation for A4 becomes 

kz[X5 + kiezX7){xs + kigXr) 



Hence A4(t) is uniquely determined as a function of x{-) and x{tf), for the right hand side of ([T9|) 
depends only on x{-),x{tf) and not on other components of A. Moreover, that X4 = implies 
jointly with the evolution equation for X4 that U2 = 0. Hence, from the condition H = const we 
deduce that A7 is constant. 
Let now 5.] = 0. Then from 

; , 3X3 + kigX7 
M = -k:,XQXj — B. 

X4 + kzXy 
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one deduces that A4 is constant, and it is given by 



62At - 
bip, — aiv 



Again, the condition H = const implies that Ay is constant, whence the whole covector is constant. 



4 Oxigen as a control 

In this section we assume (X3, X4, X5, Xg, Xy) is given such that there exists (xi, X2, X3, X4, X5, Xg, Xy) G 
T> for which a sohition of the optimal control problem in the previous section is a trajectory 
corresponding to a piecewise constant control {ul{t) , U2{t)) . 

We assume that that the oxygen affects the system and that we control oxygen concentration. 
Let U3 = O. We model /x, Vg as 



jj, = /i(x3, X4, X5, X7)a{u3) 



■De(x3,X4,X5,X7) 



1 



where a{u3) is an afHne function of Vo = ^ ^^^^ and fi, are positive functions. We assume that 
cr : R+ ^ R satisfies 

(7(0) = l>0, lim a{z) =L>1, (20) 

z— >+oo 

and cr\z) > for every z >0. We study a new optimal control problem 

max x^itf) 

0<M3<1 



(21) 



subject to 



'xs = Ful - Vg{x3,X5,X7)xe 

±4 = - Vz{x3,X4,X5,X7)xe 

X5 = Ve{x3,X4,X5,X7)x6-^^ 
X6 = il{x3, X4, X5, X7)xea{u3) 

{X7 = F, 



(22) 



and 



(X3(0), X4(0), X5(0), Xg(0), X7(0)) = {xl x^ xg, X^, X?). 



For simplicity, let us denote by x the tuple (X3, X4, X5, xe, X7) and by A the corresponding covector 
(A3, A4, A5, \q, A7). The Hamiltonian becomes 



H{x, A, ^3) = XsiFul - VgXe) + A4(Fn2 — VzXq) + xe A 



0-(w3) 



+ Xeflaius) + X7F, (23) 
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and the covector satisfies 



'As 


= Xq ( 




dVe \^ 


fe^Ae 


A4 


= X6 ( 


3^1-^4 - ^As - 






< As 


= Xq ( 


^gsF-^s + - 




3SlA6 


Ae 


= VgX 


3 + t^2A4 - t'eAs - 






A7 

V 


= Xq ( 


+ A4 - 


dVe \^ 


\ 



(24) 



The transversality conditions reads A(t/) = S/ip{x{tf)) = (0,0, 1,0,0) hence A5(t/) = 1, Xeitj) = 0. 
To maximize H{x, XjU^) with respect to u^, we compute 



dH 



xea'ius) I -As + Ae/i ) 



Clearly, XQ{t) > for every t. Notice that at the final time the maximization condition implies that 
U3{tf) maximizes XQ{t f)ve{t f) / a{u3) . If Veitf) > then u^{tf) = 0, otherwise, u^{tf) = 1. 

We say that x{-) has a singular arc on [a,h] C [0,t/], with a < 6, if = 0. As noted in the 
previous section, more degenerate situations can happen: for instance, the set of zeros of 4^ = 
can have non-isolated points without containing a whole interval. In the following theorem we 
assume that the set of zeros of |iL = is given by a whole interval and we exclude chattering 
phenomena. 

If the components AsjAg of a solution of (|24|) have only isolated zeros on [a,h\ then the only 
possibility for a singular arc, i.e., for to vanish identically, is that AsAg > and \J'^^ £ [U o'(l)] 

almost everywhere on [a, h]. In that case, ^3 = (\/^^} ■ ^^"^^ singular arc can indeed happen. 
In other words, = does not imply A5 = 0, Ae = 0. To see this, note that if there exists 
to such that A5(to) 7^ and \%{tQ) = then J^(to) 7^ 0. Similarly, if there exists to such that 
As(to) = and A6(to) / then |g(to) / 0. 'if A5(to)As(to) < then ||(to) / 0. On the 
contrary, if As(to)A5(to) > 0, then the equation ^^(to) = has a solution G [0, 1] if and only if 



A6(io)M(<o) 

If we assume that there are no chattering phenomena, then a singular arc is either of the type 
above or it satisfies As = 0, Ae = 0. In the latter case, under assumptions (pUj) on the behavior of 
a and some additional conditions on X3,X4, we show the following assertion. 

Proposition 2 Assume that the inhibition constants for glucose and xylose are different, i.e., kf^ 7^ 
k^^. Let X : [0,tf] be an optimal trajectory such that, for every t € {0,tf), x-i{t) + X4{t) > 

and both sets {t \ x^lt) = 0} and {t \ X4{t) = 0} are finite unions of intervals. If x{-) has a singular 
arc on [a, b] such that As = and Ag = 0, then A7 is constant and there are two possibilities: 

• ^3l[a,6] = 0> A3(t) = X3{a) Jq{xq{s)^^{s)) ds, A4 vanishes identically on [a,b]; 

• 3;4|[a,6] = 0, A4(t) = A4(a) Jq{x6{s)^^{s)) ds, A3 vanishes identically on [a,b]. 
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Note that we cannot prove the absence of optimal singular trajectories, not even for generic 
initial conditions. Proposition [2] only provides properties of optimal singular trajectories where 
As = and Xq = 0. 

Since Vg,Vz are given by equations ([7]), ([8]), one easily deduce that at a certain time t, Vg{t), 
respectively Vz{t), is positive if and only if x^it) > 0, respectively X4(t) > 0. Moreover, ethanol 
flux Ve is positivq^ if at least one among glucose and xylose concentrations is positive. Therefore, 
since in Section [3] we optimize the final amount of ethanol through control of feeding rates of 
xylose and glucose, it is reasonable to assume that the open loop controls ul{t),U2{t) are such that 
Vg{t) + Vz{t) > for every t. 

Proof of Proposition [2l Assume x(-) to be an optimal trajectory having a singular arc on 
[a,b] such that A5 = and Ag = 0. Since the functional to maximize ()2ip is 'ijj{x) = x^, we have 
^bi^f) = 1) which implies b < tj. Imposing conditions A5 = and Xq = 0, we obtain a linear 
homogeneous system in (A3,A4), 

1^ VgXs + VzXa = 0. 
An easy computation shows that the determinant of system (j25p is proportional to 

X3X4ikf^ - fefj, 

through a nonvanishing factor. By assumption, kf^ ~ kf^ 7^ 0. 

Case 1. Assume there exists to ^ ^] such that X3{tQ)x4,{tQ) > then (A3,A4)(to) = 0. Hence 
(-^3) ^4) I [a, 6] = 0- Imposing that A satisfies the Hamiltonian system on [a,b], since the equations 
for (A3, A4, A5, Ae) do not depend on A7 and are homogeneous in (A3, A4, A5, Ag) we deduce that 
A3,A4,A5, Ag vanish identically on the whole interval [0, ty], which contradicts A5(tj) = 1. 

Case 2. Assume now X3X4 = on [a,b]. Then, since zeros of X3,X4 are finite union of intervals 
and X3 + X4 > 0, there exists a finite decomposition [a, b] = U^|[ai, bi], {hi = Oj+i) where on each 
subinterval exactly one among X3,X4 vanishes identically and the other one is strictly positive. We 
are going to show that m is equal to 1. 

Consider a subinterval [ai,6i] where X3 = 0. Then system (|25]) implies A4 = on [aj,6j]. 
Moreover, since the evolution equation for A4 does not involve A3 and A5 = Ag = on [a, 6], we 
deduce that A4 vanishes identically on the whole interval [a,b]. Using the evolution equation for 
A3, A7 we obtain 

A3 = xe{ ^As + 7^A4 = X6^A3 

\aX3 dX3 J 0X3 

* ( dv dv \ du 

^7 = Xfs \ ^As + Tr^A4 = X67r^A3 = 0, 

\OXy OXj J OXj 

where the last identity is a consequence of Vg being proportional to X3. Hence A3|[(i- {,-] is completely 
determined as functions of x and A7 is constant on [oj , 6j] . 

^This is a consequence of the fact that growth of bacteria in the bioreactor happensx if substrates are given. Also, 
ethanol is a byproduct of the metaboHc activity of the growing bacteria. In other words, biomass growth rate /i, and 
consequently ethanol flux Ve, is positive only if sugars are given. 
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Consider now a subinterval [aj,6j] where X4 = 0. Then system (p5]) imphes A3 = on [aj,bj]. 
Using the evolution equation for A4, A7 we obtain 

M = 7^ XQAi 

0X4 

Ar = xe f I^Aa + I^A^) ^ 0, 

where the last identity is a consequence of Vz being proportional to X4. Hence, on [aj,bj] A4 is 
determined as a function of x and A3 vanishes identically. 

Let now m > 2, i.e., assume that there exists [ai, 61] C [a, b] on which X3 = and [02, 62] C [a, b] 
on which X4 = and 61 < ai. Then A4 vanishes on the whole [a,b] and A3 vanishes on [02 5^2] • 
Hence on [02,62] we have Aj = for every i = 3,4,5,6. Since the dynamics for A3,A4,A5,A6 is 
linear and homogeneous (and it does not depend on A7), this implies Aj = on the whole interval 
[0,tj]. This contradicts ^^{tf) = 1. Therefore m = 1 and the statement is proved. ■ 
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